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Abstract 



> 
ON 
\Q . 

We present the first results of an ongoing program to investigate the kinematic properties of 
high rcdshift damped Lya systems. Because damped Lya systems are widely believed to be the 
progenitors of current massive galaxies, an analysis of their kinematics allows a direct test of 
galaxy formation scenarios. Specifically, the kinematic history of protogalactic gas is a sensitive 
discriminator among competing theories of galaxy formation. 
' We use the HIRES echelle spectrograph on the Keck 10 m telescope to obtain accurate, high- 

resolution spectra of 17 damped Lya systems. We focus on unsaturated, low-ion transitions 



such as Si II 1808, since these accurately trace the velocity fields of the neutral gas dominating 
the baryonic content of the damped systems. The velocity profiles: (1) comprise multiple narrow 
components; (2) are asymmetric in that the component with strongest absorption tends to lie at 
one edge of the profile; and (3) exhibit a nearly uniform distribution of velocity widths between 
20 and 200 km s" 1 . 

In order to explain these characteristics, we consider several physical models proposed to 
explain the damped Lya phenomenon, including rapidly rotating "cold" disks, slowly rotating 
"hot" disks, massive isothermal halos, and a hydrodynamic spherical accretion model. Using 
standard Monte Carlo techniques, we run sightlines through these model systems to derive 
simulated low-ion profiles. We develop four test statistics that focus on the symmetry and 
velocity widths of the profiles to distinguish among the models. Comparing the distributions of 
test statistics from the simulated profiles with those calculated from the observed profiles, we 
determine that the models in which the damped Lya gas is distributed in galactic halos and in 
spherically infalling gas, are ruled out at more than 99.9% confidence. A model in which dwarf 
galaxies are simulated by slowly rotating "hot" disks is ruled out at 97% confidence. More 
important, we demonstrate that the Cold Dark Matter Model, as developed by Kauffmann 
(1996), is inconsistent with the damped Lya data at more than the 99.9% confidence level. This 
is because the CDM Model predicts the interception cross-section of damped Lya systems to 
be dominated by systems with rotation speeds too slow to be compatible with the data. This 
is an important result, because slow rotation speeds are generic traits of protogalaxies in most 
hierarchical cosmologies. 



1 Visiting Astronomer, W.M. Keck Telescope. The Keck Observatory is a joint facility of the University of California 
and the California Institute of Technology. 
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We find that models with disks that rotate rapidly and are thick are the only tested models 
consistent with the data at high confidence levels. A Relative Likelihood Ratio Test indicates 
disks with rotation speeds, v rot < 180 km s _1 , and scale heights, h < less than 0.1 times 
the radial scale length Rd, are ruled out at the 99% confidence level. The most likely values 
of these parameters are v ro t = 225 km s _1 and h = 0.3Rd- We also find that these disks 
must be "cold", since models in which a cc /v rot > 0.1 are ruled out with 99% confidence, where 
<r cc is the velocity dispersion of the gas. We describe an independent test of the "cold" disk 
hypothesis. The test makes use of the redshift of emission lines sometimes detected in damped 
Lya systems, as well as the absorption profiles. The test potentially distinguishes among 
damped systems that are: (1) large rotating disks detected in absorption and emission, in which 
case a systematic relation exists between emission redshift and absorption velocity profile, and 
(2) emitting galaxies surrounded by satellite galaxies detected in damped Lya absorption, in 
which case the relation between emission and absorption redshifts is random. 

Finally we emphasize a dilemma stimulated by our findings. Specifically, while the kinematics 
of the damped Lya systems strongly favor a "cold" disk-like configuration, the low metallicities 
and type II Sn abundance patterns of damped Lya systems argue for a "hot" halo-like 
configuration. We speculate on how this dilemma might be resolved. 

keywords cosmology — galaxies: evolution — galaxies: quasars — absorption lines 
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1 INTRODUCTION 



Velocity fields within protogalaxies carry important clues about the process of galaxy formation. 
The kinematic history of the protogalactic gas is of particular importance as it is a sensitive 
discriminator between competing theories of galaxy formation. In standard hierarchical models, 
low-mass subunits of dark matter and gas continuously merge to form the high-mass galaxies 
seen today [White fe Frenk 199 1| . In this scenario, the rotation speeds of galaxy disks at 
high redshifts (z > 2.5) are low compared to the rotation speeds of current massive spirals 
Kauffmann 1996]. By contrast, higher rotation speeds are predicted at large redshifts by 
models in which disks form from the coherent collapse of gas in massive dark-matter halos 
Eggen, Lynden-Bell, Sandage 19621 IFall Efstathion 1980j . 



The purpose of this paper is to test these ideas directly with data that accurately measure 
the motions of protogalactic gas. Specifically, we present an analysis of high-resolution spectra of 
metal absorption lines from a sample of damped Lya systems. The latter are high-redshift layers 
of neutral gas widely thought to be the progenitors of current galaxies W olfe 19951 IPeebles 1 993 , 



IWhit.fi Frenk 1 991 1 |Fukugita et al. 1996| . We use the HIRES Echelle spectrograph |Vogt~ 1992 



on the Keck 10 m telescope to detect the gas in absorption against the light of more distant 
QSOs. The velocity profiles of the damped systems are notable for the following properties: (1) 
they comprise multiple narrow components; (2) they are asymmetric in that the component with 
strongest absorption tends to be at one edge of the profile; and (3) the velocity intervals over which 
absorption occurs, Av, are uniformly distributed between ~ 20 and 200 km s _1 . By constructing 
statistical measures characterizing the velocity structure of the absorption profiles, we demonstrate 
that the asymmetry rules out models dominated by random motions or by spherically symmetric 
radial motions. Rather, we show the asymmetry can be caused by the rotation of "cold" gaseous 
disks with significant vertical scale heights. Furthermore, the distribution of velocity intervals 
implies high rotation speeds, v ro t ~ 225 km s _1 . 

Two types of models have been proposed to explain the damped Lya systems. The 
hierarchical schemes assume the damped systems to be slowly rotating disks embedded in low- 
mass protogalaxies |Klypin et al. 199*51 IMo fc Miralda-Escude 1994] Ma and Bertschinger 1 994 , 



IKaufimann 1 996 . Other models, that do not assume a particular cosmological context, propose 
the damped systems to be rapidly rotating disks in high-mass protogalaxies [Schiano et al. 1 990 . 
protospheroidal gas |Lu et al. 1996bj . dwarf galaxies, [York et al. 1986j . or gas falling into 
protogalaxies [Arons 1972j . We compare all the model predictions with the data using Monte Carlo 
techniques. The models are constructed by inserting gas into discrete velocity components according 
to algorithms described below. We randomly select the locations and velocities of components along 
the line of sight from parent distributions appropriate for the scenario being tested. Synthetic 
spectra are generated by computing the absorption profiles for optically thin transitions arising in 
low ions such as Si + , because the latter accurately traces the velocity structure of the neutral gas. 
We generate an ensemble of velocity profiles by assuming the model protogalaxies are randomly 
oriented in projection on the sky. We then compare the computed and observed distributions of 
statistical measures that are sensitive to the degree of asymmetry and width of the absorption 
profiles. As demonstrated below, the asymmetry in the velocity profiles rules out "hot" models 
in which a cc > v ro t at high confidence levels, where a cc is the component-component velocity 
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dispersion. This is why models with significant random motions are unlikely. We then focus 
on "cold" rotating disks in which a cc « v ro t and compare the model distribution of At; with 
the observed distribution. This particular test is crucial for evaluating model disks with different 
rotation speeds and is the test that rules out models with low v TO t- In fact, the principal conclusions 
of this paper is that low rotation speeds, characterizing the bulk of protogalactic disks, in most 
hierarchical cosmologies, are ruled out at high confidence levels. Rather, protogalactic disks appear 
to rotate rapidly. 

The paper is organized as follows. In § 2 we discuss the criteria used to select metal- line 
absorption profiles suitable for model testing. We select a sample of low-ion profiles from 17 
damped Lya systems. Optical depth spectra are constructed that will be used for tests discussed 
in the following sections. In § 3 we introduce the rotating disk model and discuss the Monte 
Carlo methods. Here we determine input parameters used to generate the simulated spectra. 
We discuss extrinsic parameters such as the number of velocity components, the velocity and 
position of each component, a cc , etc. We also discuss the selection of intrinsic parameters such 
as the internal velocity dispersion, the column density per component, etc. Together with various 
model assumptions the parameters are used to create simulated spectra with the noise and spectral 
resolution of the data. The isothermal halo, spherical infall and random models are discussed 
in § 4. The statistical tests are described in § 5 where attention is given to statistical measures 
devised to discriminate among the models. In § 6 we present the results of these tests, in particular 
demonstrating the failure of all of the models except the "cold" , rapidly rotating disk model. We 
explore the allowed parameter space for the latter model in § 7 and perform the Likelihood Ratio 
Test to determine the most likely set of values for v ro t and the thickness of the disks. An independent 
test of the rotating disk hypothesis is presented in § 8. Concluding remarks are made in § 9. 



2 PROFILE SELECTION 

This section discusses the sample of spectra used for model testing. We construct apparent optical- 
depth spectra and apply various criteria to them in order to select transitions suitable for testing. 



2.1 Data 

The data analyzed in this paper were acquired with the High Resolution Spectrograph (HIRES) 



located at the Nasmyth focus of the 10m W. M. Keck Telescope |Vogt 1992| . Table □ lists the 
coordinate name of the background QSO, observation date, exposure time, emission redshift, 
resolution and signal to noise ratio (SNR) of the damped Lya systems comprising two data sets. 
The first data set (W) comprises observations performed by one of us (AMW) as part of a long term 
project to measure the metallicity and kinematics of damped Lya systems. With the exception 
of Q1331+170, the C5 decker, that has a FWHM resolution of ~ 8 km s _1 , was used for all of 
the observations. For Q1331, we used the CI decker that gives a higher FWHM resolution of ~ 6 
km s _1 . The latter resolution and setup was used by W. L. W. Sargent who generously provided 
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Table 1: QSO AND OBSERVATIONAL DATA 
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a Data kindly provided by W. L. W. Sargent and collaborators 

b The z = 2.811 damped system toward Q0528— 2505 was omitted because its redshift exceeds that 
of the background QSO to which it may be associated. 

the second data set (S). 

All of the data were reduced with a software package kindly provided by T. Barlow. Briefly, 
the routine optimally extracts the data from its 2-D CCD format into 1-D spectra by using a bright 
standard star image as a trace. The package properly flat fields and wavelength calibrates the data 
with the quartz and ThAr calibration frames while removing cosmic rays. Finally, the data is 
continuum fitted with the IRAF package continuum. The reduction process is both accurate and 
relatively free of systematic errors. Thus the data sample is essentially homogeneous with the only 
variations arising from differences in SNR and resolution. The sample is also unbiased, because the 
objects were selected without a priori knowledge of the kinematic structure of the absorbing gas. 
Rather, the objects were selected for: (1) large H I column densit}Q N(R I) > 2xl0 20 cm" 2 , and 
(2) brightness of the background QSO, V < 19. 



lr The exceptions to this rule are the z = 1.267 system toward Q0449— 135 and the z = 1.738 system toward 
1946+770 . The low-ion column densities in both systems are large enough to indicate H I column densities, that 
have not been measured directly from Lya , exceeding the threshold iV(HI) = 2 x 10 20 cm -2 . 
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Figure 1: Velocity profiles of low-ion transitions from 17 damped Lya systems comprising our 
empirical sample. For each profile, v = km s _1 corresponds to the redshift labeled in the plot. 
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Figure 2: Average optical depth profiles (binned over 9 pixels) for the same 17 transitions presented 
in Figure 1. The solid vertical lines designate the velocity interval, Av, while the dash-dot lines 
indicate peaks and the dotted line is the 2a(f ) array. The profiles have been shifted in velocity 
space such that the left edge of the velocity interval coincides with v = and reflected if necessary 
so that the strongest peak always lies on the left edge of the interval. 
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2.2 Profile Selection Criteria 



In order to form an ensemble of velocity profiles we select the profile of a single transition from 
each damped Lya system. We first require that the line profile be free from blending with other 
absorption line profiles. This condition is met by intracomparing a set of three or more transitions 
from each damped Lya system. We next construct an apparent optical-depth spectrum, 



where I(v) and I c are the transmitted and incident intensities at velocity v. To focus on spectral 
features rather than individual pixels, we create a binned optical depth array by averaging r(v) 
over a running window nine pixels wide: 



In the last equation r(vj) is the apparent optical depth at the j th velocity pixel and the sum is over 
the four nearest neighbors of pixel i inclusive. 

In order to assure the detection of weaker features, which is crucial for establishing asymmetries 
in the line profiles, we demand that f{v p k), the f(v) of the strongest feature (i.e. the feature with 
peak optical depth), exceed 20a (f) significance. Specifically, with this criterion we can detect 
features with f(v) > 0.25xf(v p k) at 5a significance. For spectra with SNR ~ 30 this criterion 
implies the strongest feature must satisfy I(v p k)/I c < 0.6. Because saturated lines are strongly 
affected by opacity that can mask asymmetries and component structures detectable in optically 
thin transitions, we set an upper limit on f(y) by demanding I(v p k)/I c > 0.1. Therefore, we require 
0.1 < I{v p k)/Ic < 0.6 and note the upper limit will vary with SNR. Given a single transition 
satisfying this criterion, the results are independent of the ionic species or the metallicity of the 
damped systems, provided the metals trace the HI gas (i.e., the system is well mixed). 

Our analysis focuses on low-ion transitions arising in gas giving rise to damped Lya absorption 
lines. Ions such as Fe + , Si + , and Ni + are accurate tracers of neutral hydrogen when iV(H I) > 
2xl0 20 cm -2 , because the optical depth at the Lyman limit, tll > 10 3 . Moreover, the large ratios 
of A1 + /A1 ++ and Si + /Si 3+ demonstrate the gas containing low ions is mainly neutral (Prochaska 
& Wolfe 1996, 1997). The low-ion transitions selected for analysis are given in column 4 in TableEJ 
The table also gives the QSO name in column 1, the absorption redshift in column 2, the H I 
column density in column 3, the ratio f(v p k)/a(f) in column 5, and the reference for each of the 17 
systems in column 6. The corresponding intensity profiles are presented in Figure ^ The binned 
optical-depth spectra constructed from these profiles are shown in Figure HJ If necessary, we reflect 
the profiles to place the peak component at lower velocities. 

The spectra in Figures ^ and [2 are notable for several reasons. First, the unbinned spectra 
in Figure n show the profiles comprise multiple, narrow velocity components. Second, the velocity 



t(v) = Hl c /I(v)] 



(1) 




(2) 



j=i-4 
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Table 2: DAMPED LYMAN ALPHA STATISTICAL SAMPLE 
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1 j Wolfe et al.T994] 

2 jProchaska and Wolfe 1996] 

3 jProchaska and Wolfe 1997] 

4 |Ln et al. 1996b] 

5 |Lu et al. 1996a] 

6 |Lu et al. 1995a] 
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intervals over which absorption occurs are uniformly distributed between ~ 20 and 200 km s 
. Third, the velocity profiles in Figure ^ are asymmetric in that the location of the strongest 
component tends to lie at one edge of the velocity profile in 12 out of 17 cases. Fourth, the binned 
spectra in Figure|2]show evidence for a monotonic sequence of velocity components in which optical 
depth decreases with increasing velocity. To be successful, a kinematic model must quantitatively 
explain these characteristics in a natural way. 



3 ROTATING DISK MODELS 



Rotating disks have figured prominently in damped Lya research. The original survey for damped 
Lya lines was initiated in order to find rotating galactic disks at high redshifts Wol fe et al. 1 986 . 
Since then, much of the modeling has incorporated rotating disks as a working hypothesis (e.g., 
|Schiano et al. 19901 IFall Kr. Pei 19931 IKaiittma,nnT996] V The first attempt to infer the velocity 
field of the gas from the data was made by Briggs et al. (1985) who argued the red-blue asymmetry 
detected in the metal lines of a z < 1 damped absorber could arise from passage of the line of sight 
through a rotating disk. Working with higher resolution data Lanzetta & Bowen (1992) presented 
formulae for velocity profiles arising from highly inclined sightlines through a rotating spheroid, 
that mimic sightlines through a thin rotating disk, in order to model the asymmetries they found in 
three damped systems with z < 1. Using even more accurate HIRES spectra Wolfe (1995) noticed 
similar asymmetries in the metal-line profiles of five damped systems with z > 1.7. Wolfe (1995) 
assumed the velocity profiles arose in rotating disks of gas in which the gas distribution decreases 
exponentially both with radius and vertical distance from the plane. A Monte Carlo approach to 
this problem was outlined, and a preliminary description of the present work was given in Wolfe 
(1996). 

In §6 the techniques described in this section will be applied to: (1) rapidly rotating "cold" 
disks of massive galaxies; (2) slowly rotating "hot" disks of dwarf galaxies; and (3) rotating disk 
progenitors of ordinary spirals predicted by the CDM cosmology. 



3.1 Kinematic Model 



Assume the velocity profiles in damped Lya systems are explained by passage of the line of sight 
through a centrifugally supported disk in which the volume density of neutral gas is given by 



n(R, Z) = no exp 



R_ 

Rd 



\z\ 

~h 



(3) 



where R and Z are cylindrical radius and vertical displacement from midplane, uq is the central 
gas density, Rd is the radial scale length, and h is the vertical scale height. The column density of 
gas perpendicular to the disk is thus 



N ± (R) = N±(0) exp ( -R/R d 



(4) 
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where the central perpendicular column density, 



JVjl(0) = 2n h. (5) 

Because the velocity profiles in Figure ^ indicate the density field comprises a discrete rather than 
continuous distribution of gas in velocity space, n(R, Z) is proportional to the number of discrete 
clouds (i.e., velocity components) per unit volume, if, as we assume, the column density per cloud 
remains constant. Further suppose the gas rotates with a flat rotation curve given by 



v<t,(R,Z) = v rot , (6) 

where is the azimuthal component of the velocity. Although equation El breaks down at values of 
Z where the R th component of the gravitational force differs significantly from its value at Z = 0, we 
shall assume this equality is valid at all Z, and re-examine this assumption in § 9. In addition, let 
the clouds exhibit isotropic random motions characterized by a one-dimensional velocity dispersion, 
c cc . 

The profile asymmetry is the product of two effects. The "radial effect" stems from the R 
dependence of n(R, Z) along the line of sight. This produces a cloud distribution symmetric about 
the peak located where R equals the impact parameter, b, i.e., where the sightline crosses the major 
axis of the inclined disk (see Figure 01 where the major axis is denoted by the horizontal dash-dot 
line and the line of sight by the vertical dashed line). In this case the strongest absorption occurs 
at the smallest value of R where the magnitude of the rotational velocity projected along the line of 
sight is largest; i.e., the "tangent point" used in 21 cm investigations of the Galactic rotation curve 
Mih alas and Binney 198T| . As R increases, the absorption gets weaker and the magnitude of the 
projected rotational velocity decreases, thereby producing a velocity profile with strong absorption 
at one side of the profile and progressively weaker absorption toward the other. The "perpendicular 
effect" stems from the Z dependence of n(R, Z) along the line of sight. If the line of sight intersects 
midplane, Z = 0, more than a few scale heights from the major axis, the cloud distribution will 
peak at the intersection point rather than the major axis, because h < Rd- The resultant cloud 
distribution along the line of sight is lopsided with more clouds expected toward than away from 
the major axis, owing to the "radial effect". The asymmetry of this profile is the reverse of the 
previous profile, because the absorption is strongest where the magnitude of v ro t projected along 
the line of sight is smallest and weakens as the magnitude of this velocity component increases 
toward the major axis. 

Figure 01 shows examples of velocity profiles arising from these effects. The profiles designated 1 
and 2 are caused by sightlines intersecting midplane at locations 1 and 2. The profiles are produced 
by absorption in discrete clouds and are computed according to the prescription given in the next 
subsection. The asymmetry of profile 2 is dominated by the radial effect because the midplane 
intersection is so close to the major axis. Profile 1 exhibits the reverse asymmetry because the 
midplane intersection is so far from the major axis. Figure |3] also shows how the velocity profiles 
are affected by changes in b. Profiles 3—5 arise at midplane intersections with increasing b which 
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Figure 3: The center circle represents an inclined disk (top view) rotating counter clockwise with 
rotation speed v ro t = 250 km s _1 . The solid dots represent the intersection points for 5 different 
sightlines with the midplane of the disk. The solid arrows indicate the direction of the rotation 
vector. The dashed vertical line represents the sightline for intersections 1 and 2. Finally, the 
dash-dot line is the kinematic major axis, i.e. the intersection of the plane of the disk and the plane 
of the sky. The sightlines are inclined by 70 deg with respect to the normal of the disk and yield 
the profiles labeled 1—5. See the text for a discussion of the profiles. 

are displaced from the major axis by the same distance along the line of sight. The absorption 
intervals, Av, decrease with increasing b, because the gradient of v ro t projected along the line of 
sight decreases with increasing R. 

3.2 Monte Carlo Methods 

To generate a statistical ensemble of velocity profiles we adopt a Monte Carlo technique in which 
sightlines penetrate the midplanes of randomly oriented disks at random locations. The clouds 
responsible for absorption are randomly selected from the n(R, Z) distribution constrained at fixed 
impact parameter, b = Rcoscf), where 4> is the azimuthal angle in the plane of the disk. More 
specifically, we use the following prescription: 

(1) Pick input parameters N±(0), N c , ai nt , h/Rj, a cc , and v ro t, where N c is the H I 
column density per cloud, and o,- mt is the velocity dispersion internal to a cloud. For all our 
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disk models we let G{ n t = 4.3 km s . The value of Ui n t approximates the internal velocity 
dispersions of Gaussian velocity components used to model velocity profiles in damped systems 
[Wolfe et al. 19941 IProchaska and Wolfe 19961 IProchaska and Wolfe 1997] . Because the velocity 
profiles do not depend on an absolute length, we express all length scales in units of Rd- The 
remaining parameters either remain fixed during the simulation of a given ensemble of velocity 
profiles, i.e, during one run, or vary from one profile to the next. 

(2) Randomly select the angle of inclination i of the disk, where i is the angle between the 
normal to the disk and the line of sight. 

(3) Randomly draw coordinates of the midplane intersection, (X,Y), from the uniform 
projection of the inclined disk onto the plane of the sky, i.e., onto the ellipse 



X 2 Y 2 

R} + R}cos*i - ' { ) 

where the X axis coincides with the major axis of the inclined disk. We set Rf to be sufficiently 
large that sightlines outside the ellipse will not encounter jV(HI) > 2 x 10 20 cm -2 . For the values of 
iVj_(0) adopted here, we find that Rf ~ 10 Rd is satisfactory. Given the midplane intersection and 
inclination angle, the line of sight is uniquely defined by setting X equal to the impact parameter, 
b. 

(4) Compute N(X,Y), the H I column density along the line of sight, from the integral 
J n(R, Z)ds where ds is the differential length element along the line of sight and R = R(s) 

X=b 

and Z = Z(s). 

(5) Choose a minimum value for n c at the iV(HI) threshold (see below) and then for columns 
above threshold let n c be proportional to N(X, Y). After n c is selected, pick clouds randomly along 
the line of sight according to the distribution n{R, Z). 

(6) Determine the ionic column density and oscillator strength, /, that gives rise to the 
synthetic profile for the adopted transition with wavelength A. We assume A = 1808.0126 A, 
the wavelength of the Si II 1808 resonance transition, 10 12,5 cm -2 for the ionic column density of 
each cloud, and normalize / such that the product of N c f\ produces a line profile that satisfies the 
line selection criteria, i.e. having the observed f(v p k)/cr(f). Therefore, the choice of N c and A are 
arbitrary. 

(7) Compute the simulated profile by (a) solving the transfer equation for the cloud ensemble, 
(b) smoothing the spectrum to the instrumental resolution of HIRES, and (c) adding Gaussian 
noise such that the SNR of the synthetic profiles corresponds to each of the 17 observed profiles 
in Figure ^ Repeat this procedure 500 times for each observed profile in order to compute 8,500 
synthetic velocity profiles for a given run. 

Let us discuss the input parameters in more detail. The profiles are obviously sensitive to 
Vrot, h/Rd, and o- cc /v rot . In the limit a cc <C v ro t, the velocity widths are mainly affected by v ro t 
and h/Rd, because Av cx v ro t and a fixed length of sightline samples an increasing fraction of 
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Figure 4: Comparison of two line-profile fits to (a) the Si II 1808 transition from the z = 1.776 
damped system toward Q1331+170 and (b) the Cr II 2056 transition from the z = 2.040 system 
toward Q0458— 020. In both figures, the dotted histogram is the observed profile, the solid lines 
represent the individual components of the VPFIT solution, and the dash-dot line is the convolved 
final solution. Note that even though profile (a) has a smaller Av than that of the Q0458 profile, 
the VPFIT solution required 2 more clouds than the lower SNR Q0458 data. 

the rotation curve as h/R^ increases. The symmetry of the profiles is affected by cr cc /v ro t and 
h/Rdi because the systematic nature of the velocity field is diluted as the contribution made by 
random motions increases and because small values of h/R^ result in narrow symmetric profiles. 
The profiles are indirectly sensitive to iVj_(0). The portion of the disk selected for the statistical 
sample lies within the elliptical contour defined by N(X,Y) > 2xl0 20 cm -2 , hence the quantity 
iVj_(0) determines the size of the "available" disk. If N±(0) is large, the "available" disk extends 
out to many radial scalelengths, R^- In that case the impact parameter to a typical sightline will 
be large compared to Rj- This produces a small gradient in v ro t projected along the line of sight, 
that results in profiles with small Av (see sightlines 3—5 in Figure OJ). Conversely when N±(0) 
is small, larger Av will be typical. We shall refer to this effect as the "column-density kinematic 
effect". In order to reproduce the observed distribution of H I column densities, we selected -/Vj_(0) 
with values between 10 20 ' 8 and 10 21,6 cm -2 . 

On the other hand the velocity profiles are much less sensitive to n c , the number of clouds 
along the line of sight, provided v rot >> a cc ; i.e., when the kinematics are dominated by systematic 
rather than random motions. The decomposition of the absorbing gas into multiple discrete 
clouds is based on our experience with fitting model profiles to the high-resolution HIRES spectra 
(see |Wolfe et al. 1994| IProchaska and Wolfe 1996] IProchaska and Wolfe 1997] ) . Using chi square 
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fitting techniques we find that ~ five Gaussian velocity components are the minimum number 
necessary for acceptable fits, with no strict upper limit. That is, we find trials with n c = 60 yield 
test statistics (§ [SJ) similar to those trials with n c ~ 5. Although the narrowest 4 profiles in the 
observed data set can be modeled with only 1 or 2 components, these features are probably blends of 
a greater number (> 5) of components superimposed in velocity space. The case for large numbers 
of components is supported by the fact that when Av is fixed, n c increases with increasing SNR. 
In Figure |1] we show our model fits to (a) the Si II 1808 profile for the z = 1.776 damped system 
toward Q1331+170 and (b) the Cr II 2056 profile for the z = 2.040 system toward Q0458-020. 
Although Av for the Q1331 profile is below the mean for the sample, our best solution requires 
n c = 7. By contrast we need only 5 components to fit the wider, but noisier Cr II profile from the 
z = 2.04 damped system toward Q0458-020. 

Therefore, in each run we adopt n c = 5 for threshold impacts, N(X,Y) = 2 x 10 20 cm -2 , and 
then let n c be proportional to higher values of N(X, Y) for impacts within the threshold ellipse. 
The crucial point is that the results are insensitive to n c for n c > 5. The minimum value for n c 
is important, because if it were significantly smaller, velocity asymmetries could arise from strong 
Poissonian fluctuations in the strengths of the velocity components rather than from gradients in 
systematic velocity fields. We also find the results discussed below to be insensitive to crj nt , provided 
it is small compared to v ro t or a cc . 



4 MODELS WITHOUT ROTATING DISKS 

We now discuss the kinematic predictions of other models suggested to explain the damped Lya 
systems. In particular we consider models in which rotation does not dominate the velocity field. 
The Monte Carlo methods used for these tests are essentially those presented for the disk model in 



4.1 Isothermal Halos of Normal Galaxies 

The low metallicities measured in damped Lya systems, [Zn/H] ~ —1.2 |Pettini et al. 1994] . 
and the pattern of relative metal abundances in these objects Lu, Sargent, fe Barlow 1995, 



ILu et al. 1996bj bring to mind Galactic halo stars; i.e., stars that were metal enriched only by 
type II supernovae. For these reasons Lu et al. (1996b) suggest the damped systems comprise 
neutral clouds in galactic halos. The abundance patterns of damped Lya systems argue against 
metal-enriched disks of current spiral galaxies, because metal enrichment by type I as well as type 
II supernovae are required to explain their abundance patterns. We wish to test the halo hypothesis 
by modeling the cloud kinematics with an isothermal halo. 

To determine the velocity profiles we consider an isothermal halo model in which n(r), the 
gas density at 3D radius r, is proportional to the total mass density p(r) given by Table 4-1 in 
Binney and Tremaine (1987) for self gravitating isothermal spheres. The proportionality constant 
is fixed by setting the universal baryon to dark matter ratio to 0.1 |Navarro, Frenk &: White~1 995 



As a result, all the baryonic mass in galaxies is initially placed in the halo. In this model 
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Figure 5: Average optical depth profiles derived from the isothermal halo model with varying cloud 
number n c and impact parameter b. Note the profiles tend to peak near v = particularly as n c 
increases. The profiles also have very large velocity intervals, exceeding 200 km s^ 1 in the majority 
of cases. 



the mass of neutral gas in the halo cannot vary significantly between z ~ 3 and 1.5 in order 
to explain the relatively constant comoving density of neutral gas observed in this redshift 
interval :St orrie-Lombardi and W olfe 1997 . Because halos are dominated by random rather than 
systematic motions such as rotation, the asymmetries predicted for the "cold" rotating disks should 
be suppressed in this model, provided significant velocity gaps between the clouds do not exist. 
This is illustrated in Figure [5] showing eight f(v) profiles generated by a model of this type. 



4.2 Spherical Accretion onto Galaxy Bulges 

The infall of gas onto pre-existing mass concentrations is a plausible scenario for galaxy formation, 
one that may give rise to the velocity profiles in damped Lya systems. Such mass concentrations 
may have recently been discovered in the form of compact galaxy bulges with redshifts of z ~ 3.5 
Stci del et al. 1996] . The regularity of the brightness distributions further suggests these objects 
formed at even higher redshifts, and reinforces the idea that the remaining parts of the galaxy may 
have formed by infall of neutral gas onto a pre-existing bulge component [Peebles 1993j . Gaseous 
infall was first suggested by Arons (1972) who proposed that QSO absorption lines originated in 
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Figure 6: Average optical depth profiles derived from the spherical accretion model with varying 
cloud number n c , impact parameter b, and infall velocity Vi n f. Unlike all of the other models 
considered, these profiles exhibit peaks at both edges of the velocity interval. 

gas partaking in the gravitational collapse of low mass protogalaxies. Because he focused on the 
highly ionized gas responsible for the low column-density Lya forest systems, his model is not 
directly applicable to the neutral gas comprising the damped Lya systems. 

We model the infall by assuming neutral clouds responsible for the damped Lya systems are 
embedded in a hydro dynamical infall of fluid onto a spherical mass distribution. By assuming 
spherical symmetry we can adopt the analysis of Bondi (e.g. Shapiro Hi Teukolsky 1 983 ) for 



hydrodynamical spherical accretion. In this case the radial velocity and volume density of the fluid 
are given by 



i 

r s \2 



v(r) = -vo ( f ) [r « r s ] , (8) 



and 



n(r) = n — [r < r s ] , (9) 
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where vq an d no are the radial velocity and density at r = r s . The quantity r s is the transonic 
radius and is related to the ambient thermal energy per unit mass (i.e., the sound speed squared), 
ale, by 



GM 



(10) 



In this model the gas has zero net angular momentum; future work will address a scenario with 
finite net angular momentum. 

Figure El shows examples of velocity profiles produced by a typical accretion model. Notice the 
double absorption peaks located at both edges of the profile. These are caused by the maxima in 
the line-of-sight component of v(r) at displacements ± b from the point of closest approach. By 
contrast the profiles produced by "cold" rotating disks tend to produce single peaks at the profile 
edges (see Figure EJ) . 



4.3 Random Motion Models 



In order to generate a control sample, that is, a sample of velocity profiles created by velocity fields 
without inherent systematics, we investigate a model in which the velocities are randomly drawn 
from a uniform velocity distribution covering the velocity interval v = [0,v max ]. The model has a 
generic problem matching the observations in that it yields nearly the same Av for each simulated 
line profile. This is in obvious conflict with the observed profiles exhibiting a wide range of Av. 
Because the neutral volume density, n, is not described by this model, it does not yield values 
of N(H I). Therefore, for each run we arbitrarily select the number of individual clouds from the 
interval n c = [5,25], resulting in variations of iV(H I). 



5 STATISTICAL TESTS 



Figure [7| is a plot of characteristic line profiles for each of the models for a smooth distribution of 
gas (n c S> 1). Rapidly rotating disks are characterized by asymmetric distributions in optical depth 
that tend to peak at one edge of the velocity interval. The spherical accretion model produces peaks 
at both edges, the random motions model produces a flat profile without peaks, while central peaks 
or null peaks characterize "hot" models where a cc > v ro t- In this section we design test statistics 
to discriminate among these models. In particular, the tests focus on the width of the features 
in velocity space (the velocity interval), the location of the peak component within the velocity 
interval, and the degree of asymmetry evident in the profiles. 



5.1 Velocity Interval Test 



The velocity interval, Av, of an absorption profile is defined as the width of the profile in velocity 
space. A simple method for determining Av is to identify the 3a components at the lowest and 
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Figure 7: Smooth absorption profiles for the (a) Thick Rotating Disk (TRD) (b) Spherical 
Accretion (SPH) (c) Random Motion (RAN) and (d) Dwarf Galaxy (DG). Profiles derived by 
setting n c = 1000. 

highest velocities of the f(v) array and define Av as the difference in velocity between these 
components. This approach, however, gives weak outlying components the same weight as strong 
components in the main absorption complex. As a result a single, weak component displaced 200 
km s _1 from the centroid of a strong absorption complex spanning a velocity interval of 40 km s _1 
would yield Av ~ 220 km s _1 even though most of the absorption spans a much smaller Av. 

For this reason we define Av as the velocity interval weighted by optical depth. Specifically, 
define Av by first computing the integrated optical depth 



not = J r{v)dv . (11) 

We then step inward from each edge of the profile until reaching the pixel where 5% of Ttot is 
removed, and set Av equal to the velocity difference between the pixels. This definition removes 
outlying weak components that may be statistical fluctuations that are not good tracers of the 
underlying velocity field. While the routine may underestimate the true velocity interval, this is 
preferable to including misleading statistical artifacts. In addition, numerical tests (see Figure |SJ) 
reveal the removal algorithm is not sensitive to the optical depth of the absorption feature. Although 
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Figure 8: Comparison of the effects of saturation on the definition of the velocity interval by (a) 
a 3a limit in the f(v) array and (b) designating 90% of the total optical depth. Note that the 
velocity interval for the optical depth method remains constant with increasing line strength, while 
the velocity interval derived with the 3a limit clearly increases (~ 20% from the first to last frames). 



the choice of removing 10% of Ttot is somewhat arbitrary, the simulations indicate the other 
statistical tests are largely independent of the removed fraction of Ttot- 



5.2 Mean- Median Test 

One of the principal characteristics of the velocity profiles is their asymmetry, and so we devise 
a number of test statistics to discriminate asymmetric profiles from symmetric ones. The 
Mean— Median Test gauges asymmetry by measuring the relative difference between the mean 
and median velocities of the optical depth profile. The test first establishes the velocity interval 
as discussed above. The mean velocity is defined as the velocity of the midpoint of the interval 
where the edge corresponding to the smallest velocity is normalized to v = (i.e. v mean = lAv). 
The routine calculates Ttot within the velocity interval, and then determines the median velocity, 
Vmedian, as the velocity dividing Ttot in half. The formal statistic is expressed as: 



p \Vmean ^median] /"lo\ 

Jmm ~ (A^2) " 1 j 
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For a perfectly symmetric distribution of gas, f mm = 0, while more asymmetric profiles will tend 
to fmm = 1. The test is both physically motivated and relatively simple. Note that because the 
Velocity Interval routine removes 10% of Ttot symmetrically, it does not affect the Mean— Median 
Test. 



5.3 Edge-Leading Test 

The test statistic, fedg, measures the relative velocity difference between v p k, the velocity of the 
feature with peak f(v), and v me an such that 



u mean \ „\ 

Iedg = (Av/2) ■ 1 j 

Thus, f e dg = for peaks at the center of the profile, while very "edge-leading" profiles yield f e dg ~ 1- 

5.4 Two Peak Test 

A third test statistic, designed to accentuate the location of the second strongest peak in the velocity 
profiles, is given by 



\v2pk_ 

(Av/2) 



2pk = ± TaTTTa '> ( 14 ) 



where V2 P k is the velocity of the second strongest peak. The plus sign holds if the second peak 
is between the velocity of the first peak, v p k, and v m ean, otherwise the negative sign holds. 
Consequently, profiles with peaks at both edges, as would arise in the case of spherical accretion, 
or profiles in which V2 P k < v p k, that might be drawn from Gaussian velocity distributions, can 
be distinguished from the peaks arising from rotating exponential disks. The disk model tends 
to produce peaks with v p k < V2 P k < v m ean because t(v) systematically decreases with increasing 
v — Vpk- In fact, the Two— Peak Test was introduced to differentiate between the rotating disk and 
spherical infall models. In the random models the two strongest peaks should be uncorrelated and 
will occur at velocities uniformly distributed within Av. While uncorrelated, the two strongest 
peaks will be more centrally located in models for which a cc > v ro t, since the random part of the 
velocities are drawn from a Gaussian rather than a uniform velocity distribution. 

As in the Edge-Leading Test, this test first finds the strongest feature in the velocity interval. It 
then searches across Av for a second significant feature satisfying two criteria: (1) f(v2 P k) > \f(v p k)i 
and (2) f(v2 P k) is a 3er feature above an artificial continuum level set by the lowest f(vi) value 
reached between v p k and V2 P k- The first criterion prevents small, physically irrelevant, features 
from biasing the test. The second criterion is necessary to distinguish a physically distinct second 
peak from a small rise on a single feature (see Figure EJ). After the two peaks are determined, 
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Figure 9: Average optical-depth profiles for two systems, (a) one exhibiting two peaks and the 
other (b) only a single peak. The down arrows indicate peaks while the up arrow indicates the 
boundary between the two peaks. This boundary also acts as the artificial continuum level from 
which the second peak is considered a 3a feature. The dotted line is the la error array for f . 

the test measures the relative velocity of the second strongest peak from the center of the velocity 
interval. Note that if there is no significant second peak, the value for f e d g is adopted instead. 
Therefore, models preferentially yielding profiles with single peaks will have an f2 P k distribution 
biased to positive values. 

6 RESULTS 

Having generated ensembles of velocity profiles for each of the models in § 3 and § 4 we will use 
them to construct distributions of the test statistics discussed in § 5. We test the validity of the null 
hypothesis, i.e., that a given model is compatible with the data, by comparing model and empirical 
distributions with a two-sided Kolmogorov-Smirnov (KS) test. The KS test gives the probability 
Pks that the two distributions could have been drawn from the same parent population. Because 
it is the most conservative test of its kind, a low value for Pks is strong evidence that a given 
model is unlikely to represent the data. 

Figures HOfa — d show comparison plots of the 4 test statistics for each of the 6 kinematic 
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(a) Velocity Test 



(b) Mean— Median Test 




Figure 10: Plot of the Av, / e d g , fmm, and /2 P k distributions for the data and the six test models. 
The Pks values represent the probability the test distribution and the data could have been drawn 
from the same parent population. Note that every model except the TRD Model is inconsistent 
with the data for one or more statistical tests. 



Table 3: MODEL PARAMETERS 



MODEL 


N±(0) 


v rot or v 


0"cc 


n c 


h 




( cm -2 ) 


(km s _1 ) 


(km s _1 ) 


(2 x 10 20 ) 


(ro) 


TRD a 


21.2 


250 


10 


5 


0.4 


DG fe 


21.2 


50 


5-55 


5 


0.3 


CDM C 


21.2 


30-200 


10 


5 


0.3 


m d 




40 


190 


10 




SPH e 


19.4 


20 


10 


5 




RAN^ 




100 
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a Cold, Rapidly Rotating Disk Model 
b Dwarf Galaxy Model 
c Cold Dark Matter Model 
d Isothermal Halo Model 
e Spherical Accretion Model 
* Random Motions Model 



Table 4: P KS VALUES 



MODEL 


VEL 


MM 


EDG 


TRD 


TRD 


0.84 


0.37 


0.99 


0.40 


DG 


0.90 


0.034 


0.11 


0.31 


CDM 


2.0 x 10~ 4 


0.0062 


0.070 


0.49 


IH 


4.5 x 10~ 16 


0.25 


0.58 


0.0016 


SPH 


0.030 


0.90 


0.0015 


1.0 x 10" 


RAN 


0.0034 


4.5 x 10~ 4 


0.86 


4.0 x 10- 



models. After an exhaustive search of the parameter space permitted by each model, we chose 
those parameters yielding the closest agreement to the data. Table El gives the parameters used in 
the simulations, while the Prs values are listed in Table |IJ The Figure demonstrates that apart 
from "cold" disks with high rotation speeds, the remaining models are ruled out at high levels of 
confidence by at least 1 statistical test. Therefore, of all the tested models, only "cold" disks with 
rapid rotation are consistent with the damped Lya kinematics. In the following subsections, we 
will briefly discuss why the other models fail specific tests and why improved results are unlikely. 



6.1 Models with Rotating Disks 

Rotating disks are essential elements in models with "cold" disks of giant spirals, "hot" disks of 
dwarf galaxies, and the progenitors of spiral galaxies in CDM cosmologies. The models considered 
in this subsection have the properties discussed in § 3, and differ only in the choice of the parameters 
Vrot, Vcc, h/R d , n c and iVi(0). 
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6.1.1 "Cold" Disks of Giant Spirals 



In this model, velocity profiles are generated by line-of-sight impacts through a set of intrinsically 
identical exponential disks that are "cold" and rapidly rotating. Therefore, all the input parameters 
are the same for each of the 8,500 impacts comprising a given run. In order for the H I disks to 
be "cold" we let a cc = 10 km s _1 which equals the velocity dispersion measured across the H 
I disks of spiral galaxies |van der Kruit fc Shosta k 1984j . Note the results are insensitive to the 
value of a cc provided a cc < 25 km s _1 ; i.e., the disks remain "cold". For larger values of a cc , the 
model fails the Two-Peak Test and also cannot reproduce the low Av tail of the Velocity Interval 
Test distribution. For each run we set v rot = 250 km s _1 and h/Rd = 0.4. This rotation speed 
is a representative value for current giant spirals [Rubin et al. 198*5] . The scale-height is chosen to 
reproduce the observed distribution of Av. 

As this model assumes identical properties for all intercepting galaxies, it is somewhat naive. 
A more realistic model would include a distribution of rotation speeds, thicknesses and central 
column densities. However, it is valuable for testing the null hypothesis that giant spirals explain 
the kinematic state of the damped Lya systems. The figures show that the null hypothesis cannot 
be ruled out for this model, because it produces acceptable values of Pks for each of the 4 statistical 
tests. 

6.1.2 "Hot" Disks of Dwarf Galaxies 

In this model we assume the velocity profiles are generated by the disks of dwarf galaxies in which 
a cc is assumed to vary from object to object. York et al. (1986) first suggested dwarf galaxies were 
responsible for QSO absorption systems exhibiting metal lines. More recently, Pettini et al. (1994) 
and Lu et al. (1996b) advocated dwarfs as the site of damped Lya absorption. These authors were 
struck by the similarity between the low metallicities and element abundance patterns in dwarfs 
and damped Lya systems. 

The velocity profiles are generated by the same algorithms used to model the "cold" disks 
discussed in the previous subsection. That is, except for <r cc , all the input parameters are fixed for 
the profiles computed in each of 8,500 galaxy impacts. However, for these disks we assume v ro t = 50 
km s _1 which is a good representation for the rotation speeds of dwarf irregulars |Freeman"T 993 . 
We draw a cc from a distribution spanning velocities between 5 and 55 km s _1 in such a way that 
the resultant Av matches the data. The range in a cc simulates the episodic turbulence originating 
from the high frequency of supernova remnants that might occur at large redshifts. The idea is to 
let random motions in gas accelerated by supernova activity produce the observed Av. This was 
suggested earlier as the explanation behind the large Av detected in metal line absorbers as well 
as damped Lya systems [York et al. 1986| IPettini et al. 1994] , Of course, if we were to require 
a cc <C 55 km s _1 always (as observed in local gas-rich dwarfs), the model fails a priori as it cannot 
reproduce the velocity widths, Av > 100 km s" 1 . Finally, to further enhance Av we set h/Rd = 
0.3. 

Because the kinematics of the "hot" dwarf galaxy model are influenced by velocities randomly 
drawn from a Gaussian distribution, the simulated profiles tend to peak symmetrically about v = 
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km s _1 . As a result the model profiles yield discrepant / mm statistics; Pk s{f mm) = 0.037. 
Furthermore, the model is only marginally consistent with the Edge-Leading Test as Pxsifedg) 
equals 0.11. Note that the model survives the Two-Peak Test, because a large fraction of the 
profiles are produced with a single peak, which as noted above, biases the distribution to positive 
values in the f2 P k distribution. 

We wish to point out that these results are significantly dependent on both the assumed cloud 
number n c at threshold iV(HI) (taken to equal 5) and a cc . Specifically, increasing n c while holding 
a cc constant will increase the range of Av values due to the greater sampling of the random Gaussian 
distribution. Similarly, increasing a cc while holding n c constant has the same effect. At the same 
time, increasing n c reduces the level of asymmetry as the random motions begin to dominate 
the rotation kinematics. We find that the threshold n c = 5 gives the best agreement to all four 
statistical tests, and that models with threshold n c > 10 are ruled out at higher confidence levels, 
because both PKs(fmm) and Pics(fedg) would be less than 0.01. 

6.1.3 Disk Progenitors in CDM Cosmology 

Kauffmann (1996) recently combined Press-Schecter theory with the Cold Dark Matter cosmology 
(CDM) to compute the hierarchical merging of protogalactic disks embedded in dark-matter halos. 
Semi-analytic techniques were used to calculate the formation of 

H I disks as gas cools and accretes in the merging halos. With the damped Lya systems in mind, 
Kauffmann (1996) computed P(v ro t), the interception probability presented by the cross-sectional 
area of the disks, as a function of v ro t for various redshifts. The radial distribution of H I column 
densities, N±(R) was also computed for different values of v ra t- 

We calculated velocity profiles for this model as follows. For each of the 8,500 impacts 
comprising a given run we drew v ro t randomly from the P(v ro t) corresponding to the redshift z = 2.5, 
in good agreement with the average redshift, z = 2.48, of the profiles in our statistical sample. We 
then fitted exponential profiles to the plots of N±(R) by assuming -/Vj_(0) = 10 2L2 cm -2 . Kauffmann 
(1996) implicitly assumed the disks to be thin, h « Rd, and "cold", a cc « v ro t- Because v ro t is 
less than 100 km s _1 in 70% of the model disks, the predicted distribution of Av will not match 
the observed distribution which extends out to 200 km s _1 . Therefore, we assumed h = 0.3Rd to 
maximize Av without significantly changing the disk-like character of the model. Otherwise, we 
followed the same procedures used to compute velocity profiles for the other disk models. 

The results shown in Figure indicate this model is highly unlikely to represent the data. 
Owing to the predominance of slowly rotating disks, the Velocity-Interval and Mean-Median Tests 
result in Pks(Av) < 0.001 and PKs(fmm) = 0.002 which are unacceptable. Had we increased 
N±(Q) to 10 22 ' 2 cm -2 , which is closer to the value Kauffmann derived, the Av distribution would 
peak at even lower velocities, owing to the "column-density kinematic" effect, and could be ruled 
out at even higher confidence levels. Consequently we conclude this CDM model is ruled out at high 
levels of confidence. Even if we were to ignore the two highest Av values, this CDM model is still 
ruled out at > 99% c.l. The implications of this important result for theories of galaxy formation 
are discussed in § 9. 



26 



6.2 Models Without Rotating Disks 



We now turn to models without rotating disks. 

6.2.1 Isothermal Halos 

In order to mimic the kinematics of the Galaxy halo, we let the halo gas rotate with v ro t(r) = 50 
km s _1 , and assume a cc = 190 km s -1 (see Majewski 1993). As in § 13.21 we set Oi n t = 4.3 km s _1 
. For the isothermal sphere we choose the central mass density, po, in such a way that the mass 
M(r=20 kpc) = 2xlO n M , which is the mass inferred from a constant disk rotation speed v ro t = 
220 km s^ 1 out to ~ 20 kpc. The length scale of the density distribution, i.e., the King radius, 
is then determined by the relation, tq = yQo%] (AirGpo) and equals 18 kpc. However, because 
the kinematics are dominated by the random motions of the clouds, the results are essentially 
independent of po> r O) an d the assumed density profile. 

To generate simulated spectra we follow the prescription of § 13.21 except that the (X, Y) 
coordinates are drawn from a projected circle rather than a projected ellipse. The circle has a 
radius r = Rf, where Rf has the same meaning as in § 13.21 We also assume the angular momentum 
of the gas defines a preferred plane whose rotation axis is a random variable. 

Figure ITUk shows that the value of cr cc is sufficiently large to throw most of the distribution of 
Av off scale, peaking at Av ~ 600 km s _1 . As a result the model fails the Av test at a very high 
level of confidence. Furthermore, PRs(f2pk) < 0.002 owing to the large fraction of /2 P k < resulting 
from the random locations of the absorption peaks. We tried to decrease Av by reducing n c from 
10 to 5, but this increased the fraction of profiles with / 2p k < 0, and consequently produced an even 
lower PKs(f2pk) value. Because a cc » v ro t, this model, like the dwarf galaxy model, generates 
profiles exhibiting a high degree of symmetry in the limit n c ^> 1. However, in this simulation a cc 
is so large that for n c = 10, it is rare for two clouds to overlap in velocity space. Therefore, the 
strongest absorption peak occurs randomly due to the noise in the spectrum. Because this leads to 
an / e( ig distribution similar to that of the Random Motions Model, PRs(fedg) = 0.58. By the same 
reasoning, the asymmetry is accentuated such that PKs(fmm) = 0.25. However, we find that for 
n c > 15, all Isothermal Halo Models are ruled out by the Mean-Median and Edge-Leading Tests. 

6.2.2 Spherical Accretion 

We modeled the accretion flow with a central bulge having mass M ~ 10 10 Mq and gas with sound 
speed vo ph Oqo = 20 km s _1 implying r s « 110 kpc. The sound speed is physically reasonable and 
is chosen to roughly reproduce the observed Av distribution. However, it corresponds to thermal 
temperatures high enough that the gas will be collisionally ionized. Thus the gas creating damped 
Lya absorption must be discrete neutral clouds embedded in a smooth, ionized medium. 

By contrast with the other models in this section, the kinematics of the spherical infall model 
are dominated by systematic motions, i.e. radial infall, rather than random motions. However, the 
profiles are still symmetric with the peaks tending to lie at both edges of the velocity interval. As a 
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result, the Two-Peak Test distribution is dominated by negative f2 P k for all values of n c , vo and no, 
in contradiction to the empirical f2 P k distribution which is why Pr s(f 2pk) < 0.001. Furthermore, 
while the model predicts a peak in the f e a g distribution near f e d g ~ 0.75, the empirical distribution 
is more uniform between and 1.0. Thus the profiles predicted by the spherical infall model are 
actually too "edge-leading", and as a result PRs(fedg) = 0.001. Note, in order to reproduce the 
observed Av distribution we developed ad hoc distributions of vo. The combined weight of the tests 
rule out this model at very high confidence levels. 

6.2.3 Random Motion Models 

The Random Motion Model assumes a uniform deviate yielding less symmetric, more edge- leading 
profiles than those models with a Gaussian deviate. Owing to the uniform deviate, the strongest 
peak is distributed uniformly within the velocity interval. As a result the model yields acceptable 
values of PKs(fedg) f° r a wide range of cloud numbers. However, the profiles are still too symmetric 
for consistency with the Mean-Median Test for all but the smallest values of n c . And the model 
fails the Two-Peak Test and, the Velocity Interval Test for all values of n c . Therefore, this model is 
highly unlikely to represent the data. As in our treatment of the dwarf galaxies, we also considered 
a model (not plotted in Figure ITU)) in which Av is drawn uniformly from a distribution designed to 
give an acceptable distribution of Av. Although this model yields Pks(Av) ~ 1, it is ruled out at 
high confidence levels by both the Mean— Median and Two— Peak Tests. 

7 PROPERTIES OF ROTATING DISKS RESPONSIBLE FOR 
THE DAMPED Lya SYSTEMS 

Apart from rapidly rotating "cold" disks and slowly rotating "hot" disks, every model we tested is 
ruled out at higher than 99.9 % confidence by one or more of the tests. Because the "hot disks", 
i.e., the turbulent dwarf galaxies, are ruled out at 97% confidence, the only considered models 
compatible with the kinematic data are those with rapidly rotating "cold" disks. In order to place 
more rigorous bounds on parameters such as rotation speed and disk thickness, we investigate the 
parameter space of these models. In particular we perform Monte Carlo runs for a wide range of the 
parameters v ro t, h/Rd, and iVj_(0). Although we compute all four test statistics for each run, only 
the Av distributions are particularly sensitive to variations of v ro t, h/Rd, and iVj_(0). Therefore, 
we confine the discussion to the Av distributions, and to models with identical disks. 

The results are shown in Figure fTTI (left side) as iso-probability contours in the h/Rd, v ro t 
plane for N±(0) = (a) 10 20 - 8 cm~ 2 , (b) 10 212 cm~ 2 , and (c) 10 2L6 cm -2 . The contours correspond 
to selected values of Pks(Av), the KS probability discussed in §0 Because the KS test is better 
suited for ruling out rather than establishing null hypotheses, these contours indicate the range of 
v ro t and h/Rd excluded by the data; in particular, thin disks with slow rotation are unlikely. 
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Figure 11: Left Panels show isoprobability contours from the Pks(Av) tests in the v ro t, h/Rd 
plane for the cold rotating disk model with N±(0) = (a) 10 m8 cm~ 2 , (b) 10 2L2 cm -2 , and (c) 
10 21,6 cm -2 . Contour levels are drawn at Pks(Av) = 0.01, 0.05, and 0.32. Right panels show 
isoprobability contours derived from the Likelihood Ratio Test with the same contour levels. Note 
that disks with rotation speeds v ro t < 180 km s _1 or thickness log[h/Rd] < —1.0 are extremely 
unlikely. 



To find the most likely range of parameters, we turn to the Likelihood Ratio Test. Define the 
likelihood, 



L = J]KA^) , (15) 

i=l 

where p(Av{) is the probability a model yields the velocity interval, Avi, assigned to the i th member 
of a sample comprising rid damped Lya systems. We compute L max , the maximum value of L with 
respect to the parameters, i.e., we maximize the function L = L(v ro t,h/Rd, N±(0)) in the v rot , 
h/Rd plane or the v ro t, N±(0) plane. For the v ra t, h/Rd plane the maximal points are indicated 
by the x's on the right side of Figure HT1 We then compute the probability — 2 \n{L / L max ) which 
is asymptotically distributed as Ax 2 for 2 degrees of freedom in either plane, where A% 2 = \ 2 
— Xmin' an d Xmin ^ s X 1 a * ^ e maximal point. The contour levels on the right side of Figure fTTI 
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Figure 12: Isoprobability contours of Pks(Av) values in the v rot , iV(HI) plane with h/Rd = 0.3. 
The contour levels are drawn at Pks = 0.01, 0.05, and 0.32. The figure demonstrates that higher 
rotation speeds are required as N±(0) increases. 



are the confidence levels corresponding to Ax 2 > —2ln(L/L max ) for iVj_(0) = (a) 10 2a8 cm -2 , (b) 
10 21 - 2 cm" 2 , and (c) 10 21 ' 6 cm -2 . When N±(0) = 10 21 ' 6 cm -2 , L max lies on the edge of the parameter 
space we investigated and as a result determination of the "true" maximum is uncertain. We do 
not consider configurations with \og[h/Rd] > —0.25 as they do not qualify as disks. Therefore, it is 
possible that the contours with this iVj_(0) will encompass an even smaller area in this parameter 
space. 

The Likelihood Ratio Test results in several significant trends. First, a lower limit of v ro t > 180 
km s _1 is established at the 99% confidence level, independent of N±(0) and h/Rd- Upper limits 
on v ro t are less restrictive, especially for smaller values of h/Rd and larger iVj_(0). In fact there is a 
portion of parameter space allowing v ro t > 375 km s -1 at the 95% confidence level for all 3 values 
of iVj_(0). But in every case the optimal rotation speed is given by v ro t ~ 225 km s _1 . Also note 
that the 95% and 99% contours have larger areas when N±(0) = 10 2a8 cm~ 2 , a trend caused by 
the "column-density kinematic" effect discussed in § 13.21 This is further illustrated in Figure ^] 
which is an iso-probability contour plot of v ro t vs. -/Vj_(0) for h/Rd = 0.3. The gradual rise of v ro t 
with increasing -/Vj_(0) indicates larger rotation speeds are required to compensate for the narrower 
profiles produced when iVj_(0) is large. Second, the optimal h/Rd value is at 0.32. While we find a 
lower limit of h/Rd > 0.1 at the 99% confidence level, the test places no upper limit on h/Rd- This 
is because the gradient of v ro t along the line of sight increases with h/Rd until it remains nearly 
constant at h/Rd > 0.4. 

Finally, we consider the robustness of our results, in particular the optimal rotation speed, 
v ro t = 225 km s _1 , found from maximum likelihood techniques. First, we believe the rotation 
speed must be large simply because there is no evidence that the Av with the largest values are 
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Figure 13: Plot of iV(HI) vs Av. Clear squares represent the data. Dots represent simulation results 
for three runs with h/R d = 0.3, and N±(0) = (a) 10 m8 cm- 2 , (b) 10 2L2 cm- 2 , and (c) 10 2L6 cm- 2 . 



spurious outliers. Rather, the distribution out to large Av contains a tail which is real; for example, 
6 of the 17 profiles have Av > 120 km s _1 . However, our single population model with v ro t the same 
for all disks is naive in that it cannot explain the detection of a single system with Av > 225 km s _1 
. In a more realistic multi-population model, v rot would be distributed out to 350 km s _1 in order 
to explain the largest rotation speeds measured for current Sa spirals [Rubin et al. 19 85 . Such a 
model could then explain Av < 350 km s~ x . Higher values would then have to be interpreted in 
terms of multiple disks or some other mechanism. 

While the statistical tests show the general features of the "cold" disk models are consistent 
with observation, the assumption of identical -/Vj_(0) may be incorrect. Figure lT^l is a plot of observed 
iV(H I) vs Av for N±(0) = (a) 10 20 ' 8 cm" 2 , (b) 10 212 cm" 2 , and (c) 10 21 ' 6 cm" 2 . Although the data 
exhibit a possible anti-correlation between iV(H I) and Av, all three models predict a positive 
correlation between these quantities. The correlation arises because low Av are generally produced 
at large impact parameters where iV(H I) is small and large Av are produced at small impact 
parameters where iV(H I) is large. Better agreement with the data could be achieved by dropping 
the unrealistic assumption of identical iVj_(0). Thus some contribution from disks with iVj_(0) = 
10 21 - 6 cm~ 2 would help to explain data points with iV(H I) > 10 21,2 cm -2 , while disks with iVj_(0) 
_ iq20.8 cm -2 h e ip t explain the distribution of data points at iV(H I) < 10 20 ' 8 cm~ 2 . But even 
this scheme produces too many disks with Av < 60 km s" 1 and iV(H I) < 10 20 - 8 cm . These 
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might be eliminated by invoking photoionization of the outer disks where the gas density is low 
and ionization parameter high. We are currently experimenting with these and other schemes in 
order to obtain a better match with the data. 



8 AN INDEPENDENT TEST FOR RAPIDLY ROTATING 
DISKS 

Recently Lya emission, at z = 3.15, was detected from an object displaced by 2.3 arcsec from 
a sightline intersecting a damped Lya system with a similar redshift (Djorgovski et al. 1996). 
This suggests a test of the rotating disk hypothesis. Let the emission come from the center of the 
gaseous disk causing absorption. The geometry of such a configuration is illustrated in Figure 1141 
The circle represents a disk rotating with v ro t in the counterclockwise direction as indicated by 
the velocity vectors. Two possible sightline geometries are indicated by the vertical dashed lines 
that run from the observer at the bottom to the QSO at the top of the figure. The kinematic 
major axis of the disk is the dash-dot-dash line. The dots at positions 1—4 show four locations 
where the sightlines encounter the disk midplane. The metal-line velocity profiles in the figure were 
computed according to the algorithms of § 3.2 which describes why "edge-leading" profiles arising 
from midplane encounters far from the major axis, as in positions 1 and 4, exhibit asymmetries 
that mirror "edge-leading" profiles arising from midplane encounters near the major axis, as in 
positions 2 and 3. The figure shows the emission redshift, v = 0, lies at or beyond the boundaries 
of the absorption profile in every case. Furthermore, the component with peak optical depth should 
be at the profile boundary nearest to v = 0, when the point of midplane encounter is far from the 
major axis or farthest from v = 0, when the midplane encounter is nearer to the major axis. 

As a result, the disk model predicts a systematic relation between the kinematics of emission 
and absorption lines. By contrast the peak component would be randomly distributed with respect 
to v = 0, were damped Lya absorption to arise in undetected satellite galaxies bound to massive 
primary galaxies detected in emission. Suppose the speed of the satellite with respect to the primary 
is given by v sa t, and the internal kinematics of the satellite generates a profile with velocity interval, 
Av. The probability that the emission redshift of the primary is within Av is then given by 



P fl /0\ / ^ V l V sati Av<v sat , , 

Profile = (V2) | ^ ^ (16) 

Therefore, when Av « v sa t, the satellite-primary configurations are unlikely to produce emission 
redshifts within the absorption profile, and are thus indistinguishable from absorption by a single 
large disk. But when Av v sa t, the emission redshifts are likely to lie within Av. Because v sa t > 
100 km s _1 |Zaritsky White 1 995 , a significant fraction of absorption profiles with Av > 100 



km s _1 should contain the emission redshift, if the satellite/primary galaxy hypothesis is correct. 
Clearly, the combination of emission and absorption kinematics provides powerful and independent 
tests of the rotating disk model. 

Such a test was initiated by Lu et al. (1996b). They used HIRES to obtain low-ion profiles 
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Figure 14: Same as Figure 3 except here we show 2 sightlines (dashed vertical lines). In the profiles, 
the emission redshift, v = 0, lies at or beyond the boundaries of the absorption profile in every case 
and the strongest feature lies at one of the edges. 

for the z a b s = 3.15 damped Lya absorber toward Q2233+131, the same system Djorgovski et al. 
(1996) found to be associated with an emitting galaxy at the absorption redshift. Comparison 
between emission and absorption profiles shows the centroid of the emission line is located outside 
the absorption profile characterized by Av = 200 km s _1 . This is consistent with the kinematics of 
a large, rotating disk centered on the Lya emitter. Because the peak component is at the profile 
boundary farthest from v = 0, the point of midplane encounter should be near the major axis. Lu 
etal (1996b) used their measurement of At; to set a lower limit of 200 km s _1 on v ro t which was 
combined with the impact parameter, b = 30h^Q kpc, to deduce a total mass for the galaxy of M ga i 
> 3 x 10 11 /i 50 1 M Q (note this mass is somewhat higher than the |Lu et al. 1996b] estimate, since 
we assume £Iq = 0.1 rather than f^o = 1-0). On the other hand the combination of large values 
for b and Av raises potential problems for the big disk interpretation, because the column-density 
velocity effect (§ 3.2) indicates small Av when b/Rd is large. Perhaps b/Rd is small in this case 
because Rd exceeds 10 kpc, as expected from the high incidence of damped Lya systems along the 
line of sight. Moreover a large inclination angle, i > 80°, would also increase Av. 
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9 CONCLUDING REMARKS 



Our principal result is the CDM scenario for the formation of galactic disks predicts rotation speeds 
at z ~ 2.5 that are ruled out at high confidence levels by the neutral-gas kinematics in damped 
Lya systems. Specifically, the semi-analytic model of Kauffmann (1996) predicts most precursors of 
current massive spirals to be disks embedded in dark-matter halos with masses less than ~ 10 10 Mq 
and rotation speeds less than 100 km s _1 . By contrast, the likelihood ratio test in § 8 requires 
minimum rotation speeds, v ro t > 180 km s _1 , at the 99 % confidence level. The most likely value, 
v ro t = 225 km s _1 . While estimates of the masses responsible for rotation are not generally available, 
the 3 x lO 11 /i^ o 1 M0 lower limit set for the galaxy at z = 3.15 is an important constraint, because 
galaxies with such large masses rarely contain disks with significant gaseous cross sections at large 
redshifts in hierarchical scenarios j White fc Frenk 19911 iKauffmann 1996] . Rather our impression 
is that disks embedded in low-mass halos causing low rotation speeds dominate the damped Lya 
interception probabilities in most hierarchical models. We will check this impression by using our 
Monte Carlo techniques to test the predictions of other independent CDM models (see the numerical 
models of |Katz et al. 1996] ) and Hot + CDM models jKlypin et al. 19951 . Furthermore, we will 
test the recent suggestion that damped Lya absorption arises in multiple low-mass galaxies rather 
than in a single object [Rauch et al. 19961 IGardner et al. 1996j . However, until such schemes are 
shown to be consistent with the kinematic data, we shall adopt the idea of large, rapidly rotating 
disks that formed prior to z ~ 3 as a working hypothesis for the damped Lya systems. We note 
the early formation of large disks in massive galaxies is a signature of the isocurvature models 
[Peebles 1993| IPee bles 199£|, and possibly hierarchical models where fio<l and A = 0. 

Let us take a closer look at the rapidly rotating disk model. First, we assumed the rotation 
curve was flat for all radii. This is true in current spirals, owing to a coincidence at small radii 
of (i) the sum of the radial gravitational forces due to stellar disk and bulge, and (ii) the radial 
gravitational force due to the dark-matter halo Binney and Tremaine 1987| . However, in damped 



Lya systems the force due to the disk is small compared to the halo or bulge (the evidence for 
bulges in damped Lya systems stems from the recent detection of damped Lya absorption, at the 
emission redshift, in the spectrum of a Lyman-break bulge at z = 2.7 [Steidel 1996] ). because the 
surface densities of the gaseous damped Lya disks are low compared to the surface densities of 
stellar disks in current spirals. For example the central surface densities of disk-population stars 
in spirals like the Galaxy are typically 500 to 1000 M & pc~ 2 , whereas the central surface densities 
of gas in damped Lya systems are at least 10 times smaller (see [Wolfe et al. 199 5 j). As a result 
rotation curves in damped systems should behave as follows. In the absence of a bulge, the rotation 
speed due to a spherical isothermal halo rises out to the core radius and remains flat at larger radii. 
The presence of a bulge in the damped Lya system results in a steeper rise from R = because of 
the smaller core radius of the bulge. At larger radii the rotation speed will drop and then rise again 
at radii sufficiently large for the halo to dominate. Sightlines penetrating these types of disks may 
produce Av smaller than predicted by our model. We are currently constructing self-consistent 
models including bulge, disk, and halo components. Second, we assumed v ro t was independent of 
Z. In a self-consistent model the azimuthal velocity, v^(R,Z), would decrease with increasing Z 
because the R th component of the radial force decreases with Z for spherically symmetric mass 
distributions. We carried out some preliminary tests of this effect by imposing a negative Z gradient 
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on v ro t- Specifically, we let v ro t vary linearly from 225 km s -1 at Z = to 150 km s _1 at 1 kpc 
and found the gradient had only minor effects on the test statistics. Third, we are investigating the 
limits ionizing radiation places on the sizes of the gaseous disks and its effect on the test statistics. 

Having explored the gas we turn to the stars by speculating about possible stellar byproducts 
of the damped Lya disks. The kinematic evidence points to disks that: (1) are rapidly rotating, 
v ro t ~ 225 km s _1 ; (2) are "cold", a cc « v ro t; and (hi) may have large scale heights, since h > 
0.1R ra( i- By contrast the abundance measurements bring to mind "hot" configurations in which 
a cc » v ro t- This is because the low metallicity, [Zn/H] = — 1.2 |Pettini et al. 1994] or [Fe/H] = — 
1.5 |Lu et al. 1996bj . and the relative abundance patterns, suggesting gas enrichment only by type 
II supernovae |Lu et al. 1996bj . are signatures of halo population II stars. Therefore the kinematic 
and chemical evidence are in conflict. 

The conflict might be resolved if star formation in damped Lya disks results in a thick, 
metal-poor disk not unlike the thick stellar disk of the Galaxy |Gilmore, Wyse Kuijken 1 989 



The latter is characterized by scale height h ~ 1.5 kpc. The scale height is compatible with h > 
O.IRd inferred for the damped systems because Rd for current spirals ranges between 3 and 6 kpc 
[Kent 1986] . and Rd for high-z damped systems must be more than twice these values to explain 
their incidence along the line of sight [Lanzetta 1993] . However, it is not obvious that the assumed 
velocity dispersion, a cc = 10 km s _1 , can maintain matter, gas or stars, at scale heights as large 
as 1.5 kpc. At the solar neighborhood, where the surface density £ = 70 M pc -2 , stars with a cc 
~ 25 km s -1 have a scale height h = 0.4 kpc Mihalas and Binney 1981 . Because the disk models 



fail the Velocity Interval and Two Peak tests when a cc > 25 km s , large physical scale heights 
could be a serious problem for the disk model. However, the average surface density for damped 
Lya systems is only ~ IOMqpc -2 , and since a cc oc (/iS) 1 / 2 for self-gravitating disks, h ~ 1.5 kpc 
is in principal compatible with a cc < 20 km s -1 . On the other hand, the damped Lya disks may 
be self-gravitating only at R < Rd- This is because the vertical structure of the disk is controlled 
by dark matter in the halo rather than self gravity at R > Rd for exponential disks with central 
surface densities low compared to that of the Galaxy (Gunn 1982). In the dark- matter limit h 
~ (T cc R/v ro t and a cc ~ 10 km s _1 . As a result thick disks with low velocity dispersions may be 
physically plausible at any radius. While we have not run simulations in which h oc R, we see no 
reason why "flaring" disks would dilute the asymmetries required to explain the line profiles. 

Other parameters characterizing the thick stellar disk are the mass fraction, rotation speed, 
and metallicity. The mass fraction, M t hick disk/within disk ~ 0.1. This could be explained by a 
10 % decrease of £l g (z) between z = 3 and 1.5, which is consistent with the observed variation 
of ilg(z) |Storrie-Lombardi and Wolfe 1997) . The rotation speed is given by v ro t ~ 150 km s _1 . 
Although this is lower than the optimal rotation speed, v ro t = 225 km s _1 , a vertical velocity 
gradient going from v ro t = 225 km s _1 at midplane to 150 km s _1 at Z ~ h would explain the 
kinematics and results in simulated test-statistics compatible with the data, as discussed above. 
The most difficult problem has to do with metallicity, because the average [Fe/H] for the damped 
systems is between a factor of 4 and 10 below the mean metallicity of the thick disk, ([Fe/H]) = 
—0.6 |Gilmore, Wyse &: Kuijken 1989| . However, the thick disk abundance data refer only to stars 
near the solar circle, whereas the relevant metallicity is an average over the entire thick disk, which 
is unknown. As a result [Fe/H] averaged over thick stellar disks would range between —1.2 and 
—1.6 in this scenario. 
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Finally we cannot rule out the possibility that a significant fraction of the neutral gas 
experiences negligible star formation between z = 3 and 1.5. In that case the descendents of 
damped Lya systems would be extended metal-poor configurations of gas surrounding the visible 
parts of galaxies. The detection of low metallicity damped Lya systems at low redshifts, z < 1, 
provides evidence for such gas |Meyer et al. 1995 . In this case the agreement between the comoving 
density of neutral gas at z ~ 3 and the density of visible matter in current galaxies would refer to 
the inner disks where star formation occurs, presumably at z < 1.5. 
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